Multiomics reveals persistence of obesity-associated adipose tissue immune cell phenotypes during weight loss and subsequent weight regain
Figure 6
Matthew A. Cottam
Setup
library(Seurat)
library(RColorBrewer)
library(ggplot2)
library(dplyr)
library(ggpubr)
library(SingleCellExperiment)
library(DESeq2)
library(Matrix.utils)
library(magrittr)
library(stringr)
library(pheatmap)
library(tibble)
library(purrr)
library(tidyr)
library(DEGreport)
library(ggrepel)
library(reticulate)
library(viridis)
#Data file
data.integrated <- readRDS("../Inputs/IntegratedData.rds")
#Color panels
maccols <- brewer.pal(n=8, name="Blues")[c(-1,-3,-5,-7)]
monocols <- c("#ff8ade","#e324ad")
dccols <- brewer.pal(n=9, name="Greens")[-1]
tcols <- brewer.pal(n=8, name="Reds")[-1]
nkcols <- c("#876149","#6e3f22")
bcols <- brewer.pal(n=4, name="Purples")[-1]
othcols <- c("#71a157","#00c5cc","#a18e25","#b282b3")
strcols <- brewer.pal(n=4, name="Oranges")[-1]
wccols = c("#878787", "#518db6","#94cc73","#e96b53")
cols <- c(maccols,monocols,dccols,tcols,nkcols,bcols,othcols,strcols)
llcols <- c("#4292C6","#ff8ade","#238B45","#EF3B2C","#876149","#9E9AC8","#71a157","#00c5cc","#a18e25","#b282b3","#FD8D3C")
#Calculate cellspm
cellspm <- table(data.integrated$individual_mice, data.integrated$highlevel2)
propcellspm <- prop.table(cellspm, margin=1)
#Convert to a dataframe
cellspm <- as.data.frame(cellspm)
propcellspm <- as.data.frame(propcellspm)
#Change column names
colnames(cellspm) <- c("Mouse","Cluster","Counts")
colnames(propcellspm) <- c("Mouse","Cluster","Frequency")
#Add a column indicating the groups to each table
cellspm$Group <- gsub("_Hashtag.","", cellspm$Mouse)
propcellspm$Group <- gsub("_Hashtag.","", propcellspm$Mouse)
#Normalize "Counts" to "Counts per hundred cells sequenced"
cphc<-list()
cellspm <- cellspm %>% arrange(.,Mouse)
for (i in levels(cellspm$Mouse)){
x <- sum(subset(cellspm, Mouse==i)$Counts)/100
y <- subset(cellspm, Mouse==i)$Counts/x
cphc <- append(cphc, y)
}
cellspm$CountsPerHundred <- as.numeric(cphc)
#Function to calculate standard error for each given variable (cell cluster)
data_summary <- function(data, varname, groupnames){
require(plyr)
summary_func <- function(x, col){
c(mean = mean(x[[col]], na.rm=TRUE),
se = sd(x[[col]], na.rm=TRUE) / sqrt(sum(!is.na(x[[col]]))))
}
data_sum<-ddply(data, groupnames, .fun=summary_func,
varname)
names(data_sum)[names(data_sum) == 'mean'] <- varname
return(data_sum)
}
#Summarize the data
clustercounts <- data_summary(cellspm, varname="CountsPerHundred", groupnames=c("Group","Cluster"))
clusterprops <- data_summary(propcellspm, varname="Frequency", groupnames=c("Group","Cluster"))Figure 6
Figure 6A
Fig6cols <- c("#DEEBF7","#9ECAE1","#4292C6","#084594",rep("#e6e6e6",29))
DimPlot(data.integrated, split.by="orig.ident",group.by="highlevel2") +
NoLegend() +
labs(title="") &
scale_color_manual(values=Fig6cols) &
NoAxes()
Figure 6B
Fig6Macs <- subset(data.integrated, highlevel2=="Tissue Resident Macrophages" | highlevel2=="LAMs" | highlevel2=="Cycling Macrophages" |
highlevel2=="Efferocytes")Idents(Fig6Macs) <- Fig6Macs$highlevel2
DefaultAssay(Fig6Macs) <- "RNA"
VlnPlot(Fig6Macs,
features=c("Lyz2","Cst3","Adgre1","Cd68","Lgals3","Fcgr1","Itgam", "adt_MAC2","adt_CD64","adt_CD11b"),
stack=T,
flip=T,
fill.by="ident",
cols=c("#DEEBF7","#9ECAE1","#4292C6","#084594")) +
labs(x="") +
scale_x_discrete(labels=c("Tissue Resident","LAMs","Cycling","Efferocytes")) +
NoLegend() +
theme(aspect.ratio = 0.25)
Figure 6C
Fig6Macs <- subset(data.integrated, highlevel2=="Tissue Resident Macrophages" | highlevel2=="LAMs" | highlevel2=="Cycling Macrophages" |
highlevel2=="Efferocytes")Idents(Fig6Macs) <- Fig6Macs$highlevel2
DefaultAssay(Fig6Macs) <- "RNA"
VlnPlot(Fig6Macs,
features=c("Klf4","Cbr2","Stab1","Trem2","Cd9","Rgs1","Stmn1","Pclaf","Saa3","Slpi"),
stack=T,
flip=T,
fill.by="ident",
cols=c("#DEEBF7","#9ECAE1","#4292C6","#084594")) +
labs(x="") +
scale_x_discrete(labels=c("Tissue Resident","LAMs","Cycling","Efferocytes")) +
NoLegend() +
theme(aspect.ratio = 0.25) 
Figure 6D
Idents(data.integrated) <- data.integrated$orig.ident
maccellspm <- subset(cellspm,
Cluster=="Tissue Resident Macrophages" |
Cluster=="LAMs" |
Cluster=="Cycling Macrophages" |
Cluster=="Efferocytes")
maccellspm$Group <- factor(maccellspm$Group, levels = c("Lean","Obese","WL","WC"))
macclustercounts <- subset(clustercounts,
Cluster=="Tissue Resident Macrophages" |
Cluster=="LAMs" |
Cluster=="Cycling Macrophages" |
Cluster=="Efferocytes")
macclustercounts$Group <- factor(macclustercounts$Group, levels = c("Lean","Obese","WL","WC"))compare_means(CountsPerHundred ~ Group, data = subset(maccellspm, Cluster=="Tissue Resident Macrophages"), method="t.test", p.adjust.method = "bonferroni")compare_means(CountsPerHundred ~ Group, data = subset(maccellspm, Cluster=="LAMs"), method="t.test", p.adjust.method = "bonferroni")compare_means(CountsPerHundred ~ Group, data = subset(maccellspm, Cluster=="Cycling Macrophages"), method="t.test", p.adjust.method = "bonferroni")compare_means(CountsPerHundred ~ Group, data = subset(maccellspm, Cluster=="Efferocytes"), method="t.test", p.adjust.method = "bonferroni")ggboxplot(d=subset(maccellspm, Cluster=="Tissue Resident Macrophages"),
x="Group",
y="CountsPerHundred",
fill="Group",
facet.by="Cluster",
add = "mean_se") +
theme_classic() +
scale_fill_manual(values=wccols) +
theme(axis.text.x=element_text(face="bold", size=10)) +
labs(x="", y="") +
theme(axis.text.x = element_blank(), axis.ticks = element_blank()) +
facet_wrap(~Cluster) +
NoLegend() +
ylim(0,40)
ggboxplot(d=subset(maccellspm, Cluster=="LAMs"),
x="Group",
y="CountsPerHundred",
fill="Group",
facet.by="Cluster",
add = "mean_se") +
theme_classic() +
scale_fill_manual(values=wccols) +
theme(axis.text.x=element_text(face="bold", size=10)) +
labs(x="", y="") +
theme(axis.text.x = element_blank(), axis.ticks = element_blank()) +
facet_wrap(~Cluster) +
NoLegend() +
ylim(0,30)
ggboxplot(d=subset(maccellspm, Cluster=="Cycling Macrophages"),
x="Group",
y="CountsPerHundred",
fill="Group",
facet.by="Cluster",
add = "mean_se") +
theme_classic() +
scale_fill_manual(values=wccols) +
theme(axis.text.x=element_text(face="bold", size=10)) +
labs(x="", y="") +
theme(axis.text.x = element_blank(), axis.ticks = element_blank()) +
facet_wrap(~Cluster) +
NoLegend() +
ylim(0,10)
ggboxplot(d=subset(maccellspm, Cluster=="Efferocytes"),
x="Group",
y="CountsPerHundred",
fill="Group",
facet.by="Cluster",
add = "mean_se") +
theme_classic() +
scale_fill_manual(values=wccols) +
theme(axis.text.x=element_text(face="bold", size=10)) +
labs(x="", y="") +
theme(axis.text.x = element_blank(), axis.ticks = element_blank()) +
facet_wrap(~Cluster) +
NoLegend() +
ylim(0,2)
Figure 6E
FeaturePlot(subset(data.integrated, highlevel2=="Tissue Resident Macrophages"),
features="Mrc1",
split.by="orig.ident",
min.cutoff="q10",
max.cutoff="q95") &
scale_color_viridis(option="D") &
theme_void() &
NoLegend() &
xlim(-5,3) &
ylim(1,9) 
FeaturePlot(subset(data.integrated, highlevel2=="Tissue Resident Macrophages"),
features="Cd163",
split.by="orig.ident",
min.cutoff="q10",
max.cutoff="q85") &
scale_color_viridis(option="D") &
theme_void() &
NoLegend() &
xlim(-5,3) &
ylim(1,9)
#For Labels
FeaturePlot(subset(data.integrated, highlevel2=="Tissue Resident Macrophages"),
features="Cd163",
min.cutoff="q10",
max.cutoff="q85") &
scale_color_viridis(option="D") &
xlim(-5,3) &
ylim(1,9)
#For Labels
FeaturePlot(subset(data.integrated, highlevel2=="Tissue Resident Macrophages"),
features="Mrc1",
min.cutoff="q10",
max.cutoff="q95") &
scale_color_viridis(option="D") &
xlim(-5,3) &
ylim(1,9)
Figure 6F
For plotting rna velocity, I have a conda env containing scvelo (v0.2.3) and scanpy (v1.7.2). Knitting this whole file takes >~8gb of RAM. I’ve set eval=FALSE in the chunk options so that I can knit this on my local machine quickly. For running this fresh, just change to eval=TRUE.
#Prep environment for reticulate (i.e. we are running Python through R)
use_condaenv("r-velocity", required = TRUE)
conda_list()
scv <- import("scvelo")
scanpy <- import("scanpy")
matplotlib <- import("matplotlib")
plt <- import("matplotlib.pyplot", as = "plt")
#Subset for Macrophage/Monocytes clusters
Idents(data.integrated) <- data.integrated$lowlevel2
Macvelo <- subset(data.integrated, idents=c("Macrophages","Monocytes"))
#Subset for obese to reduce computation time. Looks similar for WL/WC; Lean does not have many LAMs, so trajectories did change.
Idents(Macvelo) <- Macvelo$orig.ident
Macvelo_obese <- subset(Macvelo, idents="Obese")
#To remove noise in the velocity embedding, I removed cells that were not on the same cluster trajectory (i.e. 100% disconnected).
Macvelo_cells <- CellSelector(DimPlot(Macvelo_obese))
Idents(Macvelo_obese, cells=Macvelo_cells) <- "keep"
Macvelo_obese <- subset(Macvelo_obese, idents="keep")
Idents(Macvelo_obese) <- Macvelo_obese$highlevel2
#Create adata object
spliced = Macvelo_obese@assays$spliced@counts
unspliced = Macvelo_obese@assays$unspliced@counts
row.num <- which(rownames(spliced) %in% rownames(Macvelo_obese@assays$RNA@counts))
spliced <- spliced[c(row.num),]
unspliced <- unspliced[c(row.num),]
ad <- import("anndata", convert=FALSE)
orig.ident <- Macvelo_obese$orig.ident
HTO_maxID <- Macvelo_obese$HTO_maxID
lowlevel2 <- Macvelo_obese$lowlevel2
highlevel2 <- Macvelo_obese$highlevel2
clusters <- Macvelo_obese$highlevel2
dfobs <- data.frame(orig.ident, HTO_maxID, lowlevel2, highlevel2, clusters)
rownames(dfobs) <- names(Macvelo_obese$orig.barcodes)
genes_attributes <- rownames(Macvelo_obese@assays$RNA@counts)
dfvar <- data.frame(genes_attributes)
rownames(dfvar) <- rownames(Macvelo_obese@assays$RNA@counts)
emb <- Embeddings(Macvelo_obese, "umap")
adata_Mac <- ad$AnnData(
X=t(Macvelo_obese@assays$RNA@counts),
obs=dfobs,
var=dfvar,
layers=list('spliced'=t(spliced), 'unspliced'=t(unspliced)),
obsm=list('X_umap'=emb))
adata_Mac
#Run through the scvelo pipeline and generate a dynamic velocity estimate. This takes a long time.
scv$pp$filter_genes(adata_Mac)
scv$pp$moments(adata_Mac)
scv$tl$recover_dynamics(adata_Mac)
scv$tl$velocity(adata_Mac, mode='dynamical')
scv$tl$velocity_graph(adata_Mac)
#Colors for plot
Mac.colors <- c("#DEEBF7","#9ECAE1","#4292C6","#084594","#ff8ade","#e324ad")
#Save the obese macrophage adata object
# adata_Mac$write('Reticulate/mac_obese_adata.h5ad',compression='gzip')
# adata = scv$read('Reticulate/mac_adata.h5ad')#These will open a separate window because they run through reticulate.
#save as "scvelo1.png" - This is already available in the "figures" folder.
scv$pl$velocity_embedding_stream(adata_Mac,
basis='umap',
size=100,
legend_fontsize=14,
palette=Mac.colors,
min_mass=0,
smooth=TRUE,
alpha=0.5)knitr::include_graphics("figures/scvelo1.png")
Figure 6G
In the original Figure 6 Panel G, we show only Tissue Resident Macrophages and LAMs. The code chunks below reproduce these, but also provide analysis of Efferocytes and Cycling Macrophages. In addition, we provide a code chunk for all macrophages grouped together with an additional modification to overlay gene expression on top of the MacSpectrum visualization.
Since we reproduce some of the work done for MacSpectrum and adapt it to our needs, we ask that you reach out to Dr. Beiyan Zhou and cite her lab’s work if you would like to use the code in this section.
### MacSpectrum source code was provided by Beiyan Zhou:
# Chuan Li, Antoine Menoret, Cullen Farragher, Zhengqing Ouyang, Christopher Bonin, Paul Holvoet, Anthony T. Vella, Beiyan Zhou.  Single cell transcriptomics-based MacSpectrum reveals novel macrophage activation signatures in diseases. JCI insight.  April 16 2019. PMID: 30990466.
Idents(data.integrated) <- data.integrated$orig.ident
LAMs_macspec <- subset(data.integrated, highlevel2=="LAMs")
TRMs_macspec <- subset(data.integrated, highlevel2=="Tissue Resident Macrophages")
Effs_macspec <- subset(data.integrated, highlevel2=="Efferocytes")
Cycls_macspec <- subset(data.integrated, highlevel2=="Cycling Macrophages")
all_macspec <- subset(data.integrated, lowlevel2=="Macrophages")
LAMs_counts <- as.data.frame(as.matrix(LAMs_macspec@assays$RNA@data)) %>%
mutate_all(funs(as.numeric(.))) # Macspec needs numeric data
TRMs_counts <- as.data.frame(as.matrix(TRMs_macspec@assays$RNA@data)) %>%
mutate_all(funs(as.numeric(.)))
Effs_counts <- as.data.frame(as.matrix(Effs_macspec@assays$RNA@data)) %>%
mutate_all(funs(as.numeric(.)))
Cycls_counts <- as.data.frame(as.matrix(Cycls_macspec@assays$RNA@data)) %>%
mutate_all(funs(as.numeric(.)))
all_counts <- as.data.frame(as.matrix(all_macspec@assays$RNA@data)) %>%
mutate_all(funs(as.numeric(.)))
#MacSpectrum uses ENSEMBL IDs. We use a bioMart annotation to convert our gene names to ENSEMBL IDs.
anno <- read.delim("../Inputs/Fig6/gene_to_ENS.tsv",as.is=T)
genes <- as.data.frame(rownames(LAMs_counts))
colnames(genes) <- "Gene.name"
#We lose a few genes that are not linked properly to their ENSEMBL id, and also some that are duplicated.
genes <- left_join(genes, anno) %>% drop_na() %>% distinct(Gene.stable.ID,.keep_all = TRUE) %>% distinct(Gene.name, .keep_all=TRUE)
LAMs_counts <- LAMs_counts[genes$Gene.name,]
rownames(LAMs_counts) <- genes[match(genes$Gene.name, rownames(LAMs_counts)),]$Gene.stable.ID
TRMs_counts <- TRMs_counts[genes$Gene.name,]
rownames(TRMs_counts) <- genes[match(genes$Gene.name, rownames(TRMs_counts)),]$Gene.stable.ID
Effs_counts <- Effs_counts[genes$Gene.name,]
rownames(Effs_counts) <- genes[match(genes$Gene.name, rownames(Effs_counts)),]$Gene.stable.ID
Cycls_counts <- Cycls_counts[genes$Gene.name,]
rownames(Cycls_counts) <- genes[match(genes$Gene.name, rownames(Cycls_counts)),]$Gene.stable.ID
all_counts <- all_counts[genes$Gene.name,]
rownames(all_counts) <- genes[match(genes$Gene.name, rownames(all_counts)),]$Gene.stable.ID
# Check to confirm that names of samples line up with their metadata
all(names(LAMs_macspec$individual_mice) == colnames(LAMs_counts)) ## [1] TRUEall(names(TRMs_macspec$individual_mice) == colnames(TRMs_counts)) ## [1] TRUEall(names(Effs_macspec$individual_mice) == colnames(Effs_macspec)) ## [1] TRUEall(names(Cycls_macspec$individual_mice) == colnames(Cycls_counts)) ## [1] TRUEall(names(all_macspec$individual_mice) == colnames(all_counts)) ## [1] TRUE# Store metadata
LAM_groups <- unname(LAMs_macspec$individual_mice)
TRM_groups <- unname(TRMs_macspec$individual_mice)
Eff_groups <- unname(Effs_macspec$individual_mice)
Cycl_groups <- unname(Cycls_macspec$individual_mice)
all_groups <- unname(all_macspec$individual_mice)
# Read in the index information
M1_mean <- read.table("../Inputs/Fig6/BMDM_m1_mean.txt",sep=",",header=TRUE)
M2_mean <- read.table("../Inputs/Fig6/BMDM_m2_mean.txt",sep=",",header=TRUE)
M0_mean <- read.table("../Inputs/Fig6/BMDM_m0_mean.txt",sep=",",header=TRUE)
rownames(M1_mean) <- M1_mean$GeneID
rownames(M2_mean) <- M2_mean$GeneID
rownames(M0_mean) <- M0_mean$GeneID
#
# # Save inputs - counts as a csv (instead of sparse matrix) are large.
# write.csv(LAMs_counts, file="MacSpectrum/input/LAM_genes-to-cells.csv")
# write.csv(TRMs_counts, file="MacSpectrum/input/TissueResMac_genes-to-cells.csv")
# write.csv(LAM_groups, file="MacSpectrum/input/LAM_group.csv")
# write.csv(TRM_groups, file="MacSpectrum/input/TissueResMac_group.csv")inFile <- LAMs_counts
inFile_feature <- LAM_groups
inFile <- inFile-rowMeans(inFile)
MPI_genes <- intersect(M1_mean$GeneID,rownames(inFile))
M1_mean <- M1_mean[MPI_genes,]
M2_mean <- M2_mean[MPI_genes,]
inFile_bak <- inFile
inFile <- inFile[MPI_genes,]
AMDI_genes <- intersect(M0_mean$GeneID,rownames(inFile_bak))
M0_mean <- M0_mean[AMDI_genes,]
inFile_m0 <- inFile_bak[AMDI_genes,]
#sigma of mac cells:
inFile_sigma <- 1:ncol(inFile)
total_gene_number <- nrow(inFile)
for(i in 1:ncol(inFile)) {
options(digits=9)
inFile_sigma[i] <- (sum(inFile[,i]^2)/total_gene_number)^0.5
}
#sigma of mac cells for m0:
inFile_sigma_m0 <- 1:ncol(inFile_m0)
total_gene_number <- nrow(inFile_m0)
for(i in 1:ncol(inFile_m0)) {
options(digits=9)
inFile_sigma_m0[i] <- (sum(inFile_m0[,i]^2)/total_gene_number)^0.5
}
#sigma of M0 mean:
total_gene_number <- nrow(M0_mean)
M0_sigma <- (sum(M0_mean$value^2)/total_gene_number)^0.5
#sigma of M1 mean:
total_gene_number <- nrow(M1_mean)
M1_sigma <- (sum(M1_mean$value^2)/total_gene_number)^0.5
#sigma of M2 mean:
total_gene_number <- nrow(M2_mean)
M2_sigma <- (sum(M2_mean$value^2)/total_gene_number)^0.5
#correlation of Mac to M0 mean:
total_gene_number <- nrow(inFile_m0)
inFile_Pearson_per_cell_m0 <- 1:ncol(inFile_m0)
for (j in 1:ncol(inFile_m0)){
inFile_Pearson_per_cell_m0[j] <- sum((inFile_m0[,j]/inFile_sigma_m0[j])*(M0_mean$value/M0_sigma))/total_gene_number
}
#correlation of Mac to M1 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m1 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m1[j] <- sum((inFile[,j]/inFile_sigma[j])*(M1_mean$value/M1_sigma))/total_gene_number
}
#correlation of ATM to M2 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m2 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m2[j] <- sum((inFile[,j]/inFile_sigma[j])*(M2_mean$value/M2_sigma))/total_gene_number
}
a <- 0.991414467
b <- 1
c <- -0.0185412856
x0 <- inFile_Pearson_per_cell_m1
y0 <- inFile_Pearson_per_cell_m2
d_sqr <- (a*x0+b*y0+c)^2/(a^2+b^2)
x_start <- -1
y_start <- (-a)*x_start+(-c)
x_end <- 1
y_end <- (-a)*x_end+(-c)
l <- ((x0-x_start)^2+(y0-y_start)^2-d_sqr)^0.5
l_max <- ((x_end-x_start)^2+(y_end-y_start)^2-d_sqr)^0.5
MPI <- (l-0)/(l_max-0)*100-50
AMDI <- -inFile_Pearson_per_cell_m0*50
mac_output <- data.frame(colnames(inFile),MPI,AMDI,inFile_feature,row.names=colnames(inFile))
colnames(mac_output) <- c("Samples","MPI","AMDI","Feature")
data.integrated$curr.barcodes <- rownames(data.integrated[[]])
LAMs_macspec <- left_join(x=mac_output,y=data.integrated[[c("curr.barcodes","orig.ident","HTO_maxID","lowlevel2" ,"highlevel2")]], by=c("Samples"="curr.barcodes"))
ggplot(LAMs_macspec, aes(x=MPI,y=AMDI,color=orig.ident)) +
geom_density_2d(alpha = 0.8) +
geom_point(alpha = 0.5, size=0.05) +
geom_vline(xintercept = 0) +
geom_hline(yintercept = 0) +
labs(x="MPI",y="AMDI",title="LAMs") +
theme_pubr(base_size=14) +
NoLegend() +
scale_color_manual(values=wccols) +
xlim(-16,16) +
ylim(-20,20) +
facet_wrap(~orig.ident, ncol=4) 
# This is just a figure for the annotated regions
ggplot(LAMs_macspec, aes(x=MPI,y=AMDI,color=orig.ident)) +
annotate("rect",xmin=1, xmax=16, ymin=-20,ymax=-1, fill="yellow",color=NA, alpha=0.2) +
annotate("rect",xmin=1, xmax=16, ymin=1,ymax=20, fill="red",color=NA, alpha=0.2) +
annotate("rect",xmin=-16, xmax=-1, ymin=1,ymax=20, fill="green",color=NA, alpha=0.2) +
annotate("rect",xmin=-16, xmax=-1, ymin=-20,ymax=-1, fill="blue",color=NA, alpha=0.2) +
labs(x="MPI",y="AMDI",title="") +
theme_pubr(base_size=14) +
NoLegend() +
xlim(-16,16) +
ylim(-20,20) +
geom_vline(xintercept = 0) +
geom_hline(yintercept = 0) +
theme(aspect.ratio=1)
inFile <- TRMs_counts
inFile_feature <- TRM_groups
inFile <- inFile-rowMeans(inFile)
MPI_genes <- intersect(M1_mean$GeneID,rownames(inFile))
M1_mean <- M1_mean[MPI_genes,]
M2_mean <- M2_mean[MPI_genes,]
inFile_bak <- inFile
inFile <- inFile[MPI_genes,]
AMDI_genes <- intersect(M0_mean$GeneID,rownames(inFile_bak))
M0_mean <- M0_mean[AMDI_genes,]
inFile_m0 <- inFile_bak[AMDI_genes,]
#sigma of mac cells:
inFile_sigma <- 1:ncol(inFile)
total_gene_number <- nrow(inFile)
for(i in 1:ncol(inFile)) {
options(digits=9)
inFile_sigma[i] <- (sum(inFile[,i]^2)/total_gene_number)^0.5
}
#sigma of mac cells for m0:
inFile_sigma_m0 <- 1:ncol(inFile_m0)
total_gene_number <- nrow(inFile_m0)
for(i in 1:ncol(inFile_m0)) {
options(digits=9)
inFile_sigma_m0[i] <- (sum(inFile_m0[,i]^2)/total_gene_number)^0.5
}
#sigma of M0 mean:
total_gene_number <- nrow(M0_mean)
M0_sigma <- (sum(M0_mean$value^2)/total_gene_number)^0.5
#sigma of M1 mean:
total_gene_number <- nrow(M1_mean)
M1_sigma <- (sum(M1_mean$value^2)/total_gene_number)^0.5
#sigma of M2 mean:
total_gene_number <- nrow(M2_mean)
M2_sigma <- (sum(M2_mean$value^2)/total_gene_number)^0.5
#correlation of Mac to M0 mean:
total_gene_number <- nrow(inFile_m0)
inFile_Pearson_per_cell_m0 <- 1:ncol(inFile_m0)
for (j in 1:ncol(inFile_m0)){
inFile_Pearson_per_cell_m0[j] <- sum((inFile_m0[,j]/inFile_sigma_m0[j])*(M0_mean$value/M0_sigma))/total_gene_number
}
#correlation of Mac to M1 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m1 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m1[j] <- sum((inFile[,j]/inFile_sigma[j])*(M1_mean$value/M1_sigma))/total_gene_number
}
#correlation of ATM to M2 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m2 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m2[j] <- sum((inFile[,j]/inFile_sigma[j])*(M2_mean$value/M2_sigma))/total_gene_number
}
a <- 0.991414467
b <- 1
c <- -0.0185412856
x0 <- inFile_Pearson_per_cell_m1
y0 <- inFile_Pearson_per_cell_m2
d_sqr <- (a*x0+b*y0+c)^2/(a^2+b^2)
x_start <- -1
y_start <- (-a)*x_start+(-c)
x_end <- 1
y_end <- (-a)*x_end+(-c)
l <- ((x0-x_start)^2+(y0-y_start)^2-d_sqr)^0.5
l_max <- ((x_end-x_start)^2+(y_end-y_start)^2-d_sqr)^0.5
MPI <- (l-0)/(l_max-0)*100-50
AMDI <- -inFile_Pearson_per_cell_m0*50
mac_output <- data.frame(colnames(inFile),MPI,AMDI,inFile_feature,row.names=colnames(inFile))
colnames(mac_output) <- c("Samples","MPI","AMDI","Feature")
data.integrated$curr.barcodes <- rownames(data.integrated[[]])
TRMs_macspec <- left_join(x=mac_output,y=data.integrated[[c("curr.barcodes","orig.ident","HTO_maxID","lowlevel2" ,"highlevel2")]], by=c("Samples"="curr.barcodes"))
ggplot(TRMs_macspec, aes(x=MPI,y=AMDI,color=orig.ident)) +
geom_density_2d(alpha = 0.8) +
geom_point(alpha = 0.5, size=0.05) +
geom_vline(xintercept = 0) +
geom_hline(yintercept = 0) +
labs(x="MPI",y="AMDI",title="TRMs") +
theme_pubr(base_size=14) +
NoLegend() +
scale_color_manual(values=wccols) +
xlim(-16,16) +
ylim(-20,20) +
facet_wrap(~orig.ident, ncol=4) 
inFile <- Effs_counts
inFile_feature <- Eff_groups
inFile <- inFile-rowMeans(inFile)
MPI_genes <- intersect(M1_mean$GeneID,rownames(inFile))
M1_mean <- M1_mean[MPI_genes,]
M2_mean <- M2_mean[MPI_genes,]
inFile_bak <- inFile
inFile <- inFile[MPI_genes,]
AMDI_genes <- intersect(M0_mean$GeneID,rownames(inFile_bak))
M0_mean <- M0_mean[AMDI_genes,]
inFile_m0 <- inFile_bak[AMDI_genes,]
#sigma of mac cells:
inFile_sigma <- 1:ncol(inFile)
total_gene_number <- nrow(inFile)
for(i in 1:ncol(inFile)) {
options(digits=9)
inFile_sigma[i] <- (sum(inFile[,i]^2)/total_gene_number)^0.5
}
#sigma of mac cells for m0:
inFile_sigma_m0 <- 1:ncol(inFile_m0)
total_gene_number <- nrow(inFile_m0)
for(i in 1:ncol(inFile_m0)) {
options(digits=9)
inFile_sigma_m0[i] <- (sum(inFile_m0[,i]^2)/total_gene_number)^0.5
}
#sigma of M0 mean:
total_gene_number <- nrow(M0_mean)
M0_sigma <- (sum(M0_mean$value^2)/total_gene_number)^0.5
#sigma of M1 mean:
total_gene_number <- nrow(M1_mean)
M1_sigma <- (sum(M1_mean$value^2)/total_gene_number)^0.5
#sigma of M2 mean:
total_gene_number <- nrow(M2_mean)
M2_sigma <- (sum(M2_mean$value^2)/total_gene_number)^0.5
#correlation of Mac to M0 mean:
total_gene_number <- nrow(inFile_m0)
inFile_Pearson_per_cell_m0 <- 1:ncol(inFile_m0)
for (j in 1:ncol(inFile_m0)){
inFile_Pearson_per_cell_m0[j] <- sum((inFile_m0[,j]/inFile_sigma_m0[j])*(M0_mean$value/M0_sigma))/total_gene_number
}
#correlation of Mac to M1 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m1 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m1[j] <- sum((inFile[,j]/inFile_sigma[j])*(M1_mean$value/M1_sigma))/total_gene_number
}
#correlation of ATM to M2 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m2 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m2[j] <- sum((inFile[,j]/inFile_sigma[j])*(M2_mean$value/M2_sigma))/total_gene_number
}
a <- 0.991414467
b <- 1
c <- -0.0185412856
x0 <- inFile_Pearson_per_cell_m1
y0 <- inFile_Pearson_per_cell_m2
d_sqr <- (a*x0+b*y0+c)^2/(a^2+b^2)
x_start <- -1
y_start <- (-a)*x_start+(-c)
x_end <- 1
y_end <- (-a)*x_end+(-c)
l <- ((x0-x_start)^2+(y0-y_start)^2-d_sqr)^0.5
l_max <- ((x_end-x_start)^2+(y_end-y_start)^2-d_sqr)^0.5
MPI <- (l-0)/(l_max-0)*100-50
AMDI <- -inFile_Pearson_per_cell_m0*50
mac_output <- data.frame(colnames(inFile),MPI,AMDI,inFile_feature,row.names=colnames(inFile))
colnames(mac_output) <- c("Samples","MPI","AMDI","Feature")
data.integrated$curr.barcodes <- rownames(data.integrated[[]])
Effs_macspec <- left_join(x=mac_output,y=data.integrated[[c("curr.barcodes","orig.ident","HTO_maxID","lowlevel2" ,"highlevel2")]], by=c("Samples"="curr.barcodes"))
ggplot(Effs_macspec, aes(x=MPI,y=AMDI,color=orig.ident)) +
geom_density_2d(alpha = 0.8) +
geom_point(alpha = 0.5, size=0.05) +
geom_vline(xintercept = 0) +
geom_hline(yintercept = 0) +
labs(x="MPI",y="AMDI",title="Effs") +
theme_pubr(base_size=14) +
NoLegend() +
scale_color_manual(values=wccols) +
xlim(-16,16) +
ylim(-20,20) +
facet_wrap(~orig.ident, ncol=4) 
inFile <- Cycls_counts
inFile_feature <- Cycl_groups
inFile <- inFile-rowMeans(inFile)
MPI_genes <- intersect(M1_mean$GeneID,rownames(inFile))
M1_mean <- M1_mean[MPI_genes,]
M2_mean <- M2_mean[MPI_genes,]
inFile_bak <- inFile
inFile <- inFile[MPI_genes,]
AMDI_genes <- intersect(M0_mean$GeneID,rownames(inFile_bak))
M0_mean <- M0_mean[AMDI_genes,]
inFile_m0 <- inFile_bak[AMDI_genes,]
#sigma of mac cells:
inFile_sigma <- 1:ncol(inFile)
total_gene_number <- nrow(inFile)
for(i in 1:ncol(inFile)) {
options(digits=9)
inFile_sigma[i] <- (sum(inFile[,i]^2)/total_gene_number)^0.5
}
#sigma of mac cells for m0:
inFile_sigma_m0 <- 1:ncol(inFile_m0)
total_gene_number <- nrow(inFile_m0)
for(i in 1:ncol(inFile_m0)) {
options(digits=9)
inFile_sigma_m0[i] <- (sum(inFile_m0[,i]^2)/total_gene_number)^0.5
}
#sigma of M0 mean:
total_gene_number <- nrow(M0_mean)
M0_sigma <- (sum(M0_mean$value^2)/total_gene_number)^0.5
#sigma of M1 mean:
total_gene_number <- nrow(M1_mean)
M1_sigma <- (sum(M1_mean$value^2)/total_gene_number)^0.5
#sigma of M2 mean:
total_gene_number <- nrow(M2_mean)
M2_sigma <- (sum(M2_mean$value^2)/total_gene_number)^0.5
#correlation of Mac to M0 mean:
total_gene_number <- nrow(inFile_m0)
inFile_Pearson_per_cell_m0 <- 1:ncol(inFile_m0)
for (j in 1:ncol(inFile_m0)){
inFile_Pearson_per_cell_m0[j] <- sum((inFile_m0[,j]/inFile_sigma_m0[j])*(M0_mean$value/M0_sigma))/total_gene_number
}
#correlation of Mac to M1 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m1 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m1[j] <- sum((inFile[,j]/inFile_sigma[j])*(M1_mean$value/M1_sigma))/total_gene_number
}
#correlation of ATM to M2 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m2 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m2[j] <- sum((inFile[,j]/inFile_sigma[j])*(M2_mean$value/M2_sigma))/total_gene_number
}
a <- 0.991414467
b <- 1
c <- -0.0185412856
x0 <- inFile_Pearson_per_cell_m1
y0 <- inFile_Pearson_per_cell_m2
d_sqr <- (a*x0+b*y0+c)^2/(a^2+b^2)
x_start <- -1
y_start <- (-a)*x_start+(-c)
x_end <- 1
y_end <- (-a)*x_end+(-c)
l <- ((x0-x_start)^2+(y0-y_start)^2-d_sqr)^0.5
l_max <- ((x_end-x_start)^2+(y_end-y_start)^2-d_sqr)^0.5
MPI <- (l-0)/(l_max-0)*100-50
AMDI <- -inFile_Pearson_per_cell_m0*50
mac_output <- data.frame(colnames(inFile),MPI,AMDI,inFile_feature,row.names=colnames(inFile))
colnames(mac_output) <- c("Samples","MPI","AMDI","Feature")
data.integrated$curr.barcodes <- rownames(data.integrated[[]])
Cycls_macspec <- left_join(x=mac_output,y=data.integrated[[c("curr.barcodes","orig.ident","HTO_maxID","lowlevel2" ,"highlevel2")]], by=c("Samples"="curr.barcodes"))
ggplot(Cycls_macspec, aes(x=MPI,y=AMDI,color=orig.ident)) +
geom_density_2d(alpha = 0.8) +
geom_point(alpha = 0.5, size=0.05) +
geom_vline(xintercept = 0) +
geom_hline(yintercept = 0) +
labs(x="MPI",y="AMDI",title="Cycls") +
theme_pubr(base_size=14) +
NoLegend() +
scale_color_manual(values=wccols) +
xlim(-16,16) +
ylim(-20,20) +
facet_wrap(~orig.ident, ncol=4) 
inFile <- all_counts
inFile_feature <- all_groups
inFile <- inFile-rowMeans(inFile)
MPI_genes <- intersect(M1_mean$GeneID,rownames(inFile))
M1_mean <- M1_mean[MPI_genes,]
M2_mean <- M2_mean[MPI_genes,]
inFile_bak <- inFile
inFile <- inFile[MPI_genes,]
AMDI_genes <- intersect(M0_mean$GeneID,rownames(inFile_bak))
M0_mean <- M0_mean[AMDI_genes,]
inFile_m0 <- inFile_bak[AMDI_genes,]
#sigma of mac cells:
inFile_sigma <- 1:ncol(inFile)
total_gene_number <- nrow(inFile)
for(i in 1:ncol(inFile)) {
options(digits=9)
inFile_sigma[i] <- (sum(inFile[,i]^2)/total_gene_number)^0.5
}
#sigma of mac cells for m0:
inFile_sigma_m0 <- 1:ncol(inFile_m0)
total_gene_number <- nrow(inFile_m0)
for(i in 1:ncol(inFile_m0)) {
options(digits=9)
inFile_sigma_m0[i] <- (sum(inFile_m0[,i]^2)/total_gene_number)^0.5
}
#sigma of M0 mean:
total_gene_number <- nrow(M0_mean)
M0_sigma <- (sum(M0_mean$value^2)/total_gene_number)^0.5
#sigma of M1 mean:
total_gene_number <- nrow(M1_mean)
M1_sigma <- (sum(M1_mean$value^2)/total_gene_number)^0.5
#sigma of M2 mean:
total_gene_number <- nrow(M2_mean)
M2_sigma <- (sum(M2_mean$value^2)/total_gene_number)^0.5
#correlation of Mac to M0 mean:
total_gene_number <- nrow(inFile_m0)
inFile_Pearson_per_cell_m0 <- 1:ncol(inFile_m0)
for (j in 1:ncol(inFile_m0)){
inFile_Pearson_per_cell_m0[j] <- sum((inFile_m0[,j]/inFile_sigma_m0[j])*(M0_mean$value/M0_sigma))/total_gene_number
}
#correlation of Mac to M1 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m1 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m1[j] <- sum((inFile[,j]/inFile_sigma[j])*(M1_mean$value/M1_sigma))/total_gene_number
}
#correlation of ATM to M2 mean:
total_gene_number <- nrow(M2_mean)
inFile_Pearson_per_cell_m2 <- 1:ncol(inFile)
for (j in 1:ncol(inFile)){
inFile_Pearson_per_cell_m2[j] <- sum((inFile[,j]/inFile_sigma[j])*(M2_mean$value/M2_sigma))/total_gene_number
}
a <- 0.991414467
b <- 1
c <- -0.0185412856
x0 <- inFile_Pearson_per_cell_m1
y0 <- inFile_Pearson_per_cell_m2
d_sqr <- (a*x0+b*y0+c)^2/(a^2+b^2)
x_start <- -1
y_start <- (-a)*x_start+(-c)
x_end <- 1
y_end <- (-a)*x_end+(-c)
l <- ((x0-x_start)^2+(y0-y_start)^2-d_sqr)^0.5
l_max <- ((x_end-x_start)^2+(y_end-y_start)^2-d_sqr)^0.5
MPI <- (l-0)/(l_max-0)*100-50
AMDI <- -inFile_Pearson_per_cell_m0*50
mac_output <- data.frame(colnames(inFile),MPI,AMDI,inFile_feature,row.names=colnames(inFile))
colnames(mac_output) <- c("Samples","MPI","AMDI","Feature")
data.integrated$curr.barcodes <- rownames(data.integrated[[]])
alls_macspec <- left_join(x=mac_output,y=data.integrated[[c("curr.barcodes","orig.ident","HTO_maxID","lowlevel2" ,"highlevel2")]], by=c("Samples"="curr.barcodes"))
ggplot(alls_macspec, aes(x=MPI,y=AMDI,color=orig.ident)) +
geom_density_2d(alpha = 0.8) +
geom_point(alpha = 0.5, size=0.05) +
geom_vline(xintercept = 0) +
geom_hline(yintercept = 0) +
labs(x="MPI",y="AMDI",title="alls") +
theme_pubr(base_size=14) +
NoLegend() +
scale_color_manual(values=wccols) +
xlim(-16,16) +
ylim(-20,20) +
facet_wrap(~orig.ident, ncol=4) 
#For plotting other gene information:
gene <- "Retnla"
alls_macspec$goi <- unname(all_macspec@assays$RNA@data[gene,])
ggplot(alls_macspec, aes(x=MPI,y=AMDI,color=goi)) +
geom_density_2d(alpha = 0.8) +
geom_point(alpha = 0.5, size=0.05) +
geom_vline(xintercept = 0) +
geom_hline(yintercept = 0) +
labs(x="MPI",y="AMDI",title=gene) +
theme_pubr(base_size=14) +
scale_color_viridis(option="D") +
theme(legend.position="right", aspect.ratio = 1) +
xlim(-16,16) +
ylim(-20,20) +
facet_wrap(~orig.ident, ncol=4)
gene <- "Il1b"
alls_macspec$goi <- unname(all_macspec@assays$RNA@data[gene,])
ggplot(alls_macspec, aes(x=MPI,y=AMDI,color=goi)) +
geom_density_2d(alpha = 0.8) +
geom_point(alpha = 0.5, size=0.05) +
geom_vline(xintercept = 0) +
geom_hline(yintercept = 0) +
labs(x="MPI",y="AMDI",title=gene) +
theme_pubr(base_size=14) +
scale_color_viridis(option="D") +
theme(legend.position="right", aspect.ratio = 1) +
xlim(-16,16) +
ylim(-20,20) +
facet_wrap(~orig.ident, ncol=4)
Session Info
sessionInfo()## R version 4.1.0 (2021-05-18)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur 10.16
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRblas.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] parallel stats4 stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] plyr_1.8.6 viridis_0.6.1
## [3] viridisLite_0.4.0 reticulate_1.20
## [5] ggrepel_0.9.1 DEGreport_1.28.0
## [7] tidyr_1.1.3 purrr_0.3.4
## [9] tibble_3.1.2 pheatmap_1.0.12
## [11] stringr_1.4.0 magrittr_2.0.1
## [13] Matrix.utils_0.9.8 Matrix_1.3-4
## [15] DESeq2_1.32.0 SingleCellExperiment_1.14.1
## [17] SummarizedExperiment_1.22.0 Biobase_2.52.0
## [19] GenomicRanges_1.44.0 GenomeInfoDb_1.28.0
## [21] IRanges_2.26.0 S4Vectors_0.30.0
## [23] BiocGenerics_0.38.0 MatrixGenerics_1.4.0
## [25] matrixStats_0.59.0 ggpubr_0.4.0
## [27] dplyr_1.0.7 ggplot2_3.3.5
## [29] RColorBrewer_1.1-2 SeuratObject_4.0.2
## [31] Seurat_4.0.3
##
## loaded via a namespace (and not attached):
## [1] utf8_1.2.1 tidyselect_1.1.1
## [3] RSQLite_2.2.7 AnnotationDbi_1.54.1
## [5] htmlwidgets_1.5.3 grid_4.1.0
## [7] BiocParallel_1.26.0 Rtsne_0.15
## [9] munsell_0.5.0 codetools_0.2-18
## [11] ica_1.0-2 future_1.21.0
## [13] miniUI_0.1.1.1 withr_2.4.2
## [15] colorspace_2.0-2 highr_0.9
## [17] knitr_1.33 ROCR_1.0-11
## [19] ggsignif_0.6.2 tensor_1.5
## [21] listenv_0.8.0 labeling_0.4.2
## [23] lasso2_1.2-21.1 GenomeInfoDbData_1.2.6
## [25] mnormt_2.0.2 polyclip_1.10-0
## [27] farver_2.1.0 bit64_4.0.5
## [29] parallelly_1.27.0 vctrs_0.3.8
## [31] generics_0.1.0 xfun_0.24
## [33] doParallel_1.0.16 R6_2.5.0
## [35] clue_0.3-59 isoband_0.2.5
## [37] locfit_1.5-9.4 reshape_0.8.8
## [39] bitops_1.0-7 spatstat.utils_2.2-0
## [41] cachem_1.0.5 DelayedArray_0.18.0
## [43] assertthat_0.2.1 promises_1.2.0.1
## [45] scales_1.1.1 gtable_0.3.0
## [47] Cairo_1.5-12.2 globals_0.14.0
## [49] goftest_1.2-2 rlang_0.4.11
## [51] genefilter_1.74.0 GlobalOptions_0.1.2
## [53] splines_4.1.0 rstatix_0.7.0
## [55] lazyeval_0.2.2 spatstat.geom_2.2-2
## [57] broom_0.7.7 yaml_2.2.1
## [59] reshape2_1.4.4 abind_1.4-5
## [61] backports_1.2.1 httpuv_1.6.1
## [63] tools_4.1.0 psych_2.1.6
## [65] logging_0.10-108 ellipsis_0.3.2
## [67] spatstat.core_2.3-0 jquerylib_0.1.4
## [69] ggdendro_0.1.22 ggridges_0.5.3
## [71] Rcpp_1.0.7 zlibbioc_1.38.0
## [73] RCurl_1.98-1.3 rpart_4.1-15
## [75] deldir_0.2-10 GetoptLong_1.0.5
## [77] pbapply_1.4-3 cowplot_1.1.1
## [79] zoo_1.8-9 grr_0.9.5
## [81] haven_2.4.1 cluster_2.1.2
## [83] data.table_1.14.0 scattermore_0.7
## [85] openxlsx_4.2.4 circlize_0.4.13
## [87] lmtest_0.9-38 RANN_2.6.1
## [89] tmvnsim_1.0-2 fitdistrplus_1.1-5
## [91] hms_1.1.0 patchwork_1.1.1
## [93] mime_0.11 evaluate_0.14
## [95] xtable_1.8-4 XML_3.99-0.6
## [97] rio_0.5.27 readxl_1.3.1
## [99] shape_1.4.6 gridExtra_2.3
## [101] compiler_4.1.0 KernSmooth_2.23-20
## [103] crayon_1.4.1 htmltools_0.5.1.1
## [105] mgcv_1.8-36 later_1.2.0
## [107] geneplotter_1.70.0 DBI_1.1.1
## [109] ComplexHeatmap_2.8.0 MASS_7.3-54
## [111] car_3.0-10 igraph_1.2.6
## [113] forcats_0.5.1 pkgconfig_2.0.3
## [115] foreign_0.8-81 plotly_4.9.4.1
## [117] spatstat.sparse_2.0-0 foreach_1.5.1
## [119] annotate_1.70.0 bslib_0.2.5.1
## [121] XVector_0.32.0 digest_0.6.27
## [123] ConsensusClusterPlus_1.56.0 sctransform_0.3.2
## [125] RcppAnnoy_0.0.18 spatstat.data_2.1-0
## [127] Biostrings_2.60.1 rmarkdown_2.9
## [129] cellranger_1.1.0 leiden_0.3.8
## [131] edgeR_3.34.0 uwot_0.1.10
## [133] curl_4.3.2 shiny_1.6.0
## [135] rjson_0.2.20 lifecycle_1.0.0
## [137] nlme_3.1-152 jsonlite_1.7.2
## [139] carData_3.0-4 limma_3.48.0
## [141] fansi_0.5.0 pillar_1.6.1
## [143] lattice_0.20-44 Nozzle.R1_1.1-1
## [145] KEGGREST_1.32.0 fastmap_1.1.0
## [147] httr_1.4.2 survival_3.2-11
## [149] glue_1.4.2 zip_2.2.0
## [151] iterators_1.0.13 png_0.1-7
## [153] bit_4.0.4 stringi_1.7.3
## [155] sass_0.4.0 blob_1.2.1
## [157] memoise_2.0.0 irlba_2.3.3
## [159] future.apply_1.7.0